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SUMMARY 

A two-dimensional, nonaxisymmetric Euler solution in a geometry representative of a jet engine 
configuration without blades is presented. The domain, including internal and external flow, 
is covered with a multiblock grid. In order to construct this grid, a domain decomposition 
technique is used to subdivide the domain, and smooth grids are dimensioned and placed in 
each block. The Euler solution is verified by examining five theoretical properties. The result 
demonstrates techniques for performing numerical solutions in complex geometries and provides 
a foundation for complete engine throughflow calculations. 


NOMENCLATURE 


C D 

drag coefficient 

t 

time 

C L 

lift coefficient 

U,V 

velocity components 

D() 

artificial dissipation 

vol 

volume 

E 

total energy 

w 

vector of dependent variables 

n 

surface normal 

a 

angle of attack 

H 

total enthalpy 

7 

ratio of specific heats 

P 

pressure 

P 

density 

QO 

flux integral approximation 

cn 

stage coefficients 

s 

surface area 

n 

arbitrary region 


INTRODUCTION 

Numerical simulations have provided valuable insights into a wide range of important problems 
in science and engineering. In aeronautical engineering these applications range from funda- 
mental problems in fluid mechanics to practical problems in design. Simulations have helped 
reduce design time while improving designs, and they are now an integral part of the design 
process. 

Traditionally, jet engine research and development has involved a substantial experimental com- 
ponent. Engineering practice for engine development in the 1980’s involved preliminary analysis 
and computer modeling of the first-order physical effects followed by extensive experimental 
testing of engine components in isolation. This effort was followed by engine integration based 
on component performance and culminated in ground testing of the integrated engine. The 
process is expedited by experience and extensive databases of test data, which allow estimates 
of physical effects including endwall boundary layer losses and shock losses. However, this de- 
velopment and integration can cost hundreds of millions of dollars, and it can be years before 
the viability of a design is firmly established. In order to reduce the risks of this process, de- 
sign changes are often incremental and build on experience, and because of these risks, radical 
designs are difficult to fully test. 


l 



Computational methods are available that complement experimental techniques in the analy- 
sis of engine flow phenomena. For isolated blade rows inviscid and viscous codes have been 
developed. Dawes 1 as well as Chima and Yokota 2 have reported codes for three-dimensional 
viscous flow in a blade row. Whitfield et al. 3 have demonstrated techniques for time-accurate 
simulations of counterrotating propfans. 

In multistage compressors and turbines today’s algorithms and computers dictate the use of 
physical modeling because of the geometrical complexities of relatively moving blade rows and 
the physical effects of their interactions. A foundation of compressor throughflow calculations 
is the work of Wu 4 who has presented an approach to solving quasi-three-dimensional flow by 
using a meridional surface coupled with a blade-to-blade surface. Jennions and Stowe 5 - 6 have 
demonstrated a throughflow calculation for a turbine blade. Wisler et al. 7 have mentioned the 
use of circumferentially averaged throughflow calculations for evaluating compressor designs. 
Recently, a number of techniques have been introduced for three-dimensional multistage calcu- 
lations. Denton 8 has presented an Euler calculation for multistage turbomachines, and Dawes 9 
has reported a three-dimensional viscous technique. To model flow in a multistage machine, 
Adamczyk et al. 10 have developed a combination of space- and time-averaging operators in the 
spirit of Reynolds averaging. 

Despite the complexities of simulating components, they do not operate in isolation and their 
interactions are important. Yet, in engine development, components are integrated on the 
basis of their performance in isolation. With predictive tools largely limited to one-dimensional 
analysis of complete engines, designs must be conservative even when there is a need for less 
conservative designs that require consideration of these interaction effects. 

The goal of the current work is to lay a foundation for complete engine throughflow calculations. 
As a step in this task the current calculation deals with the geometrical and numerical problem 
of solving for the external and internal flow in the inlet, bypass duct, core duct, and nozzle of the 
engine configuration. The geometry is not axisymmetric in that there is a uniform depth over 
the computational domain instead of radial variation. Given this foundation, determination of 
the axisymmetric unbladed solution will follow. 11 A further step is modeling blade effects in the 
solution. 


First, the task of determining the contours for this geometrical configuration is detailed. Then 
the techniques for finding a multiblock grid within the domain are presented. Finally, the 
numerical methods are explained and solution results and their verification are presented. 


GEOMETRY MODELING 

The geometry for this problem is based on the Energy Efficient Engine (E 3 ), which was designed 
and tested by General Electric for NASA in the early 1980’s. In particular, the geometry is based 
on the ground test of the integrated engine. The engine was designed for use in commercial 
jets in the 1990’s. 

Several modifications were made to the geometry; the most important was neglecting blade 
effects in the compressor and turbine sections. The compressor does work and compresses the 
air passing through the engine core, and for the actual engine the pressure ratio is 23:1. After 
passing through the combustor, where there is an enthalpy increase, the flow expands through 
turbine blades, which remove energy from the flow to drive the compressor and the fan. In this 


solution these effects were neglected, and the compressor, combustor, and turbine sections were 
replaced with a duct. However, a duct with a comparable area ratio through these sections will 
produce sonic flow at the throat even at relatively low Mach numbers. A shock at the throat 
limits mass flow through the engine core, and changes the streamlines throughout the engine. 
In order to reduce this problem the core area ratio was changed to 2:1 by raising the splitter 
plate. 

Several other modifications were made to the geometry contours. Since the geometry was for a 
ground test, an external nacelle was designed by using standard design practices. Furthermore, 
a small splitter island near the core bypass was excluded. Since surface pressure is very sensitive 
to the smoothness of the surface contours, every effort was made to ensure continuity of the 
surfaces and their derivatives. In particular, least-squares curve fitting to a set of sine and 
cosine basis functions has been used to remove high-frequency errors in the coordinate data. 


GRID GENERATION 

The grid on which the calculation is done is a multiblock grid that consists of block-structured 
grids which cover the domain without overlapping. The interfaces between blocks are matched 
so that coordinate lines pass from one block to another continuously without a slope discon- 
tinuity. These grids are generated by using the TOPOS program 12 , which finds a domain 
decomposition, dimensions the grids within the blocks, and places smooth grids within each 
block. 

The domain decomposition is determined with an algorithm that finds boundary-conforming 
regions within a domain. In a way analogous to how the skin of a balloon conforms to neigh- 
boring walls when blown up in a confined space, this algorithm finds topologically rectangular 
regions within the domain. If this region is removed from the domain to produce a truncated 
domain, the algorithm may be applied again to find a further block. By repeated application 
of the algorithm, the domain is reduced until it is covered. 

Once a domain decomposition has been found, transformations may be applied to cut and split 
the decomposition and change its topology. The transformations are isomorphic with respect to 
the data structures used to represent the domain decomposition. These transformations, with 
their isomorphic properties, provide a valuable tool for manipulating domain decompositions. 

The dimensions of each grid block are found so that coordinate lines continue through each block 
interface. Since the block topology is globally unstructured, this interface property places a 
nontrivial constraint on the grid dimensions. However, these properties may be formulated as 
a system of underconstrained linear equations with constant, integer coefficients that require 
positive integer solutions. If the system of equations admits a solution, this solution can be 
found with the simplex linear programming algorithm. 13 Furthermore, the simplex method 
specifies some of the parameters in this system of equations that may be varied to adapt the 
grid. Smooth grids are then placed in these dimensioned blocks by using a Coons patch followed 
by smoothing with an elliptic grid generator. 
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NUMERICAL METHODS 


The two-dimensional Euler equations model inviscid, compressible flow and are given by 


where F = { f,g ) and 
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Here 7 = 1.4. The numerical approximation and solution techniques are based on the FL052 
program of Jameson et al . 14 The discretization of these equations on the grid takes U to be a cell, 
and the flux through each face is approximated from the centroid-based dependent variables as 


I 

face 
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where M and N are cells adjacent to the face. The discretized equations for each cell are 
advanced to a steady state by a multi-stage scheme, 

w n) = Wi + ai A£[Q(u><) + D(wi)] 
w ^ — Wi + + D(w (1) )J 

w (3) = Wi + <t 3 A£[Q(w (2) ) + (3) 

= Wi + a4At{Q(w^) + D(w (1) )] 

W i+ 1 = Wi -I- a 5 At{Q(w (i) ) + Z)(w (1) )j 

where ^ are stage coefficients, Q{w) is the convective flux approximation for the cell given in 
equation (2), and D(w) is an artificial dissipation. The artificial dissipation consists of third- 
order dissipation, which stabilizes the time-stepping scheme, and first-order dissipation, which 
is switched on near shocks to capture them. 


Three types of boundaries can occur in a grid. First, at solid surfaces a free-slip, tangential 
flow condition is applied. Second, at the truncation boundary in the far field the free-stream 
conditions are matched to internal conditions by matching incoming and outgoing Riemann 
invariants. Third, at the interface boundary between blocks information must be exchanged. 
Since coordinate lines pass through these boundaries continuously, the continuation of one grid 
into its neighbor is trivial. Each grid in computational space has a layer of cells surrounding 
it, where dependent variable values from the continuation onto the neighboring grid may be 
placed. Then the convective flux Q(w) and the artificial dissipation D{w) are calculated, as usual, 
in each cell within the grid. By using this scheme and passing dependent variables across block 
interfaces at each stage of equation (3), the solution is not influenced by these block interfaces 
regardless of how the domain is decomposed. 
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Several techniques are used to enhance the convergence of equation (3). Since there are no blades 
that do work on the flow, the solution is isenthalpic and enthalpy damping 15 may be used to 
enhance convergence. Since a steady-state solution to equation (1) is sought, the maximum 
local time step permitted by the CFL condition is used. Consequently, the simulation is not 
limited by a small time step determined by the smallest cell in the grid. 

A multigrid algorithm 16 is used to accelerate convergence further. The multigrid levels are 
constructed by finding the coarsest grid that does not have excessive stretching or high aspect 
ratios due to the nonuniformity of the decomposition topology. Each dimension of this coarse 
grid is multiplied by a power of 2 to construct finer grid levels. In the current calculation the 
grid topology and coarsest grid allow three grid levels to be used and this improves convergence 
considerably. For each grid level added, the convergence rate, measured by the average density 
residual, increases by about a factor of 2. 


RESULTS 

The grid on which the calculation has been done is shown in figure 1, where the bold lines 
indicate the interfaces between blocks. The grid contains 16800 cells in 10 blocks that cover 
the external field, inlet, bypass duct, core duct, and nozzle of the nonaxisymmetric engine 
configuration. 



Figure 1: Inner grid for jet engine configuration. 
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The geometry is symmetric about the meanline of the hub, but the grid is not since there is no 
symmetry condition applied to the grid between the two halves. With a symmetric grid at zero 
angle of attack the solution would be exactly symmetric and measures of the solution, such as 
the lift coefficient would exactly cancel. However, with an asymmetric grid the solution will 
not necessarily be symmetric, and the lift coefficient will not necessarily be zero. Consequently, 
this grid asymmetry can be exploited to verify the resolution of the solution. 

The simulation was run 2500 iterations at Mach 0.3 and 0° angle of attack. The convergence 
history is shown in figure 2. The single-grid convergence history is also shown to demonstrate 
the speedup due to multigrid. 



Figure 2: Convergence histories of average density residual 
with two multigrid levels (solid line) and without multigrid (dashed line). 

There are no known experimental results for comparison, but the simulation may be verified 
on the basis of five theoretical quantities. First, the solution conserves mass and energy. Con- 
servation of mass is directly demonstrated in figure 2 by the convergence of the average density 
residual |^| over the domain. Conservation of energy follows from convergence of the average 
total energy residual |^| over the domain. Second, since the angle of attack, a, is zero and 
the geometry is symmetric, the lift coefficient, C L , should be zero. It is calculated to be C L 
= -0.0031. The deviation from zero is due to the asymmetry of the grid, the truncation error 
in the convective flux approximation Q{w), and the artificial dissipation D{w). Third, since the 
solution is subcritical, the drag coefficient C D should be zero. It is calculated to be C D = 0.0045. 
Again, the deviation from zero is due to the asymmetry of the grid, error in the approximation 
of Q(W), and the artificial dissipation D(w). Fourth, since inviscid, adiabatic flow is isentropic, 
deviations from constant entropy are a further test. Over the surfaces of the body the entropy 
deviations are less than 0.5 percent. The large deviations are at the leading edge of the splitter 
plate with its small leading-edge radius, at the curved inlet to the core passage, and at the 
leading edge of the nacelles. 
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Fifth, the asymmetry of the grid does not guarantee a symmetric solution. Figure 3 shows the 
pressure distributions over the symmetric surfaces of the components. Surface pressure is sensi- 
tive to errors, and the result shows symmetry since surface pressure overplots on the symmetric 
surfaces. The largest deviation is in the core passage and is attributed to an accumulation of 
error as the flow passes through the core duct. 




Figure 3: Normalized pressure distributions over the symmetric surfaces. Clockwise from top left, 
the components are the hub, the splitter plate, and the nacelle. 

CONCLUSIONS 

A two-dimensional, nonaxisymmetric Euler solution for a jet engine geometry has been pre- 
sented. The multiblock grid is generated by using techniques that find a domain decomposition 
and a grid within the domain. The solution has been verified on the basis of five theoretical 
quantities. These results demonstrate a foundation on which further simulations can be built, 
including solutions for axisymmetric configurations and solutions with blade effects included. 
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